
#ifndef __STROKE2D_H_
#define __STROKE2D_H_

#include "../typedefs.h"
#include <list>
#include <vector>
#include <GL/glut.h>

using namespace std;

class Stroke {
   
public:
   ///
   Stroke() {
      
   }
   
   ///
   unsigned size() {
      return positions.size();
   }
   
   ///
   void clear() {
      positions.clear();
   }
   
   ///
   void push_back(const Point2 &p) {
      positions.push_back(p);
   }
   
   void genLocalCoordinates() {
      
      Point2 o = positions.front(); // the origin
      Point2 x = positions.back();
      Vector2 u = x - o;
      u = u/sqrt(u*u);
      Vector2 v(-u[1],u[0]);
      
      Point2 y = o + 2*v;
      axisX = Segment2(o,x);
      axisY = Segment2(o,y);
      
      localpos.resize(positions.size());
      alpha.resize(positions.size(),0.0);
      beta.resize(positions.size(),0.0);
      d.resize(positions.size(),0.0);
      
      float tx, dx, ty, dy;
      unsigned int n = positions.size();
      for (unsigned i=0; i<n; ++i) {
         const Point2 &p = positions[i];
         //~ dx = computeDistanceSegment(p, axisX, tx);
         //~ dy = computeDistanceSegment(p, axisY, ty);
         //~ localpos[i] = Point2(dy,dx);
         dx = (p-o)*u;
         dy = (p-o)*v;
         localpos[i] = Point2(dx,dy);
      }
      
      beta[0] = 0.0;
      for (unsigned i=1; i<n; ++i) {
         const Point2 &pb = localpos[i-1]; // point before
         const Point2 &p  = localpos[i];
         beta[i] = (p[1]-pb[1])/(p[0]-pb[0]);
      }
      
      //~ cout<<"beta ..."<<endl;
      //~ for (uint i=0; i<beta.size(); ++i) {
         //~ cout<<i<<" "<<beta[i]<<endl;
      //~ }
      
      d[0] = 0.0;
      for (unsigned i=1; i<n; ++i) {
         const Point2 &pb = localpos[i-1]; // point before
         const Point2 &p = localpos[i];
         const Vector2 v = p - pb;
         d[i] = d[i-1] + sqrt(v*v);
      }
      
      for (unsigned i=0; i<n; ++i) {
         d[i] = d[i]/d[n-1];
      }
      
      alpha[0]   = beta[1];
      alpha[n-1] = beta[n-1];
      for (unsigned i=1; i<n-1; ++i)
         alpha[i] = (d[i]*beta[i] + d[i+1]*beta[i+1])/(d[i]+d[i+1]);
      
      for (unsigned i=0; i<n; ++i) {
         cout<<"i "<<i<<" :d: "<<d[i]<<" :alpha: "<<alpha[i]<<" :beta: "<<beta[i]<<endl;
      }
   }
   
   Number computeDistanceSegment(const Point2 &v, const Segment2 &segment, Number &t) {
      
      const Point2 &a = segment.source();
      const Point2 &b = segment.target();
      
      const Vector2 &ba = b-a;
      const Vector2 &va = v-a;
      const Vector2 &vb = v-b;
      
      const Number &norma2aibi = ba*ba;
      const Number &dotvaibiai = va*ba;
      
      Number dminvaibi;
      if (dotvaibiai<=0)                  { dminvaibi = sqrt(va*va); t=1.0;}
      else if (dotvaibiai>=norma2aibi)    { dminvaibi = sqrt(vb*vb); t=0.0;}
      else                                { 
         dminvaibi = sqrt(va*va - ((dotvaibiai*dotvaibiai)/norma2aibi)); 
         t = 1.0-dotvaibiai/norma2aibi;
      }
      
      return dminvaibi;
   }
   
   /// Simplifies p so that each edge has squared length ~= minLen 
   /// @param p boundary to be simplified
   /// @param minLen square of length of smallest admissible edge
   void subsampleStroke(double minLen=30*30) {
      
      vector<Point2>::iterator firstp = positions.begin(), endp = positions.end();
      Point2 prevpt = positions.front();
      vector<Point2> newPositions;
      
      for (vector<Point2>::iterator ip = positions.begin(); ip != positions.end(); ++ip) {
         Point2 pt = *ip;
         double sqdist = squared_distance(prevpt, pt);
         if (sqdist >= minLen) {
            supersampleEdge(newPositions, prevpt, pt, sqdist/minLen);
            newPositions.push_back(pt);
            prevpt = pt;
         }
      }
      
      Point2 pfirst = positions.front();
      Point2 pback = positions.back();
      positions = newPositions;
      if (squared_distance(positions.front(), pfirst) > 10) 
         positions.insert(positions.begin(), pfirst);
      if (squared_distance(positions.back(), pback) > 10) 
         positions.push_back(pback);
      
      genLocalCoordinates();
   }
   
   ///
   void draw() {
      
      // Draw the polygonal line
      glPushAttrib(GL_LIGHTING_BIT | GL_POINT_BIT | GL_LINE_BIT);
      glDisable(GL_DEPTH_TEST);
      glDisable(GL_LIGHTING);
      glEnable(GL_POINT_SMOOTH);
      glEnable(GL_LINE_SMOOTH);
      glEnable(GL_BLEND);
      
      glColor4f(1.0,0.5,0.0,1.0);
      glLineWidth(2.0);
      glPointSize(8.0);
      setMatrices();
      vector<Point2>::const_iterator pi = positions.begin(), pbeyond = positions.end();
      glBegin(GL_POINTS);
      while (pi != pbeyond) {
         Point2 pwindow = *pi++;
         Point3 pworld = screenToWorld(pwindow[0], pwindow[1]);
         glVertex3f(pworld[0], pworld[1], pworld[2]);
      }
      glEnd();
      
      glColor4f(0.3,0.3,0.0,1.0);
      glBegin(GL_LINES);
         Point3 o = screenToWorld(axisX.source()[0], axisX.source()[1]);
         Point3 x = screenToWorld(axisX.target()[0], axisX.target()[1]);
         Point3 y = screenToWorld(axisY.target()[0], axisY.target()[1]);
         glVertex3f(o[0], o[1], o[2]);
         glVertex3f(x[0], x[1], x[2]);
         glVertex3f(o[0], o[1], o[2]);
         glVertex3f(y[0], y[1], y[2]);
      glEnd();
      
      {
      glColor4f(0.5,0.,0.0,1.0);
      vector<Point2>::const_iterator pi = localpos.begin(), pbeyond = localpos.end();
      glBegin(GL_POINTS);
      while (pi != pbeyond) {
         Point2 pwindow = *pi++;
         Point3 pworld = screenToWorld(pwindow[0], pwindow[1]);
         glVertex3f(pworld[0], pworld[1], pworld[2]);
      }
      glEnd();
      }
   
      glPopAttrib();
   }
   
private:
   
   /// Returns the world coordinates of a point in screen space
   Point3 screenToWorld(int x, int y) {
      
      GLdouble xw, yw, zw;
      gluUnProject(x, y, 0.0, mMtx, pMtx, viewport, &xw, &yw, &zw);
      return Point3(xw, yw, zw);
   }
   
   /// 
   void setMatrices() {
      
      glGetIntegerv(GL_VIEWPORT, viewport);
      glGetDoublev(GL_PROJECTION_MATRIX, pMtx);
      glGetDoublev(GL_MODELVIEW_MATRIX, mMtx);
   }
   
   /// Supersamples an edge.
   /// @param newp new boundary where edge points will be added
   /// @param p0 first endpoint of edge
   /// @param p1 second endpoint of edge
   /// @param n  ratio between the squared length of the edge and the squared 
   /// length of each desired edge segment
   void supersampleEdge (vector<Point2>& newpos, Point2 p0, Point2 p1, double n) {
      
      if (n > 4) {
         Point2 mid = midpoint(p0, p1);
         supersampleEdge(newpos, p0, mid, n/4);
         newpos.push_back(mid);
         supersampleEdge(newpos, mid, p1, n/4);
      }
   }
   
public:
   vector<float>     d;
   vector<float>     alpha;     ///
private:   
   vector<Point2>    positions; /// positions of the stroke in window coordinates
   vector<Point2>    localpos;  /// positions of the stroke in local coordinates 
   vector<float>     beta;      /// 

   
   Segment2          axisX;
   Segment2          axisY;
   
   GLdouble          mMtx[16], pMtx[16];  /// modelview/projection matrix
   GLint             viewport[4];         /// the viewport
};

#endif
